/*
 *  Scalar.cpp
 *  T3nsors
 *
 *  Created by Michael Barriault on 09-10-07.
 *  Copyright 2010 University of Guelph. All rights reserved.
 *
 */

#include "Scalar.h"
#include "Macros.h"
#include "Others.h"
#include <cmath>
#include <cstdarg>
#include <cstdio>
#include <cstdlib>
using std::string;

Scalar::Scalar(void) {
	return;
}

Scalar::Scalar(Domain iO, std::string iid) {
	O = iO;
	id = iid;
	T.assign(O(), 0.0);
}

Scalar::~Scalar(void) {
	T.clear();
}

double& Scalar::operator[](int i) {
	if (!(T[i] == T[i])) {
		printf("\nGot a NaN here, item %d of object %s\n", i, id.c_str());
		exit(-1);
	}
	return T[i];
}

double& Scalar::operator()(int a, ...) {
	va_list va;
	va_start(va, a);
	int t=0;
	FOR(i,O.N) {
		int A = (a>=0) ? a : O[i]+a;
		for ( int j=i-1; j>=0; j-- )
			A *= O[j];
		t += A;
		a = va_arg(va, int);
	}
	va_end(va);
	return T[t];
}

Scalar operator-(Scalar x) {
	Scalar y(x.O,(string)"-"+x.id);
	FOR(o,x.O())
		y[o] = -x[o];
	return y;
}

Scalar operator+(Scalar x, Scalar y) {
	CheckDims(x.O,y.O);
	Scalar z(x.O,x.id+"+"+y.id);
	FOR(o,x.O())
		z[o] = x[o] + y[o];
	return z;
}

Scalar operator-(Scalar x, Scalar y) {
	CheckDims(x.O,y.O);
	Scalar z(x.O,x.id+"-"+y.id);
	FOR(o,x.O())
		z[o] = x[o] - y[o];
	return z;
}

Scalar operator*(Scalar x, Scalar y) {
	CheckDims(x.O,y.O);
	Scalar z(x.O,x.id+"*"+y.id);
	FOR(o,x.O())
		z[o] = x[o] * y[o];
	return z;
}

Scalar operator/(Scalar x, Scalar y) {
	CheckDims(x.O,y.O);
	Scalar z(x.O,x.id+"/"+y.id);
	FOR(o,x.O())
		z[o] = x[o] / y[o];
	return z;
}

Scalar operator^(Scalar x, Scalar y) {
	CheckDims(x.O,y.O);
	Scalar z(x.O,x.id+"^"+y.id);
	FOR(o,x.O())
		z[o] = pow(x[o], y[o]);
	return z;
}

Scalar operator+(Scalar x, double y) {
	Scalar z(x.O,x.id+"+"+ftos(y));
	FOR(o,x.O())
		z[o] = x[o] + y;
	return z;
}

Scalar operator-(Scalar x, double y) {
	Scalar z(x.O,x.id+"-"+ftos(y));
	FOR(o,x.O())
		z[o] = x[o] - y;
	return z;
}

Scalar operator*(Scalar x, double y) {
	Scalar z(x.O,x.id+"*"+ftos(y));
	FOR(o,x.O())
		z[o] = x[o] * y;
	return z;
}

Scalar operator/(Scalar x, double y) {
	Scalar z(x.O,x.id+"/"+ftos(y));
	FOR(o,x.O())
		z[o] = x[o] / y;
	return z;
}

Scalar operator^(Scalar x, double y) {
	Scalar z(x.O,x.id+"^"+ftos(y));
	FOR(o,x.O())
		z[o] = pow(x[o] , y);
	return z;
}

Scalar operator+(double x, Scalar y) {
	return y+x;
}

Scalar operator-(double x, Scalar y) {
	return -(y-x);
}

Scalar operator*(double x, Scalar y) {
	return y*x;
}

Scalar operator/(double x, Scalar y) {
	Scalar z(y.O,ftos(x)+"/"+y.id);
	FOR(o,y.O())
		z[o] = x / y[o];
	return z;
}

Scalar operator^(double x, Scalar y) {
	Scalar z(y.O,ftos(x)+"^"+y.id);
	FOR(o,y.O())
		z[o] = pow(x , y[o]);
	return z;
}

Scalar sin(Scalar x) {
	Scalar y(x.O,(string)"sin"+x.id);
	FOR(o,x.O())
		y[o] = sin(x[o]);
	return y;
}

Scalar cos(Scalar x) {
	Scalar y(x.O,(string)"cos"+x.id);
	FOR(o,x.O())
		y[o] = cos(x[o]);
	return y;
}

Scalar sin2(Scalar x) {
	Scalar y(x.O,(string)"sin2"+x.id);
	FOR(o,x.O())
		y[o] = pow(sin(x[o]),2.);
	return y;
}

Scalar sqrt(Scalar x) {
	Scalar y(x.O,(string)"sqrt"+x.id);
	FOR(o,x.O())
		y[o] = sqrt(x[o]);
	return y;
}

double l2(Scalar x, double p) {
	double l2x = 0;
	FOR(o,x.O()) l2x += pow(x[o],p);
	l2x *= pow(x.O.d(), x.O.N);
	FOR(n,x.O.N) l2x /= (x.O.max(n)-x.O.min(n));
	l2x = pow(fabs(l2x),1./p);
	return l2x;
}
